t=1;
% interest
R=@(s) cos(exp(-s));
% volatility
sig=@(s) 1+s-s;
%ic
y0=zeros(21,1);
yic=1;
y0(1)=yic;
% solver
tspan = linspace(0,t,100);
options = odeset('AbsTol',1e-12,'RelTol',1e-6);
[testt,testy] = ode45(@testrdforce2,tspan,y0,[],R,sig);